#encoding=utf8
import sys
import math
import json

def calLatLngDistance(lng1, lat1, lng2, lat2):
    rlat1 = math.radians(lat1)
    rlat2 = math.radians(lat2)
    a = rlat1 - rlat2
    b = math.radians(lng1) - math.radians(lng2)
    s = 2 * math.asin(math.sqrt(
    math.pow(math.sin(a/2), 2) + math.cos(rlat1)*math.cos(rlat2)*math.pow(math.sin(b/2), 2)))
    return s * 6372797.0

# latlngList   [ (106.45485,29.562632,-43,"wjscy") ,(106.45484,29.562631,-82,"wjscy") ,(106.45483,29.562632,-30,"wjscy") ....]
# e  邻域
# minpts 最少个数
# calDistance 距离计算公式
def dbscan(latlngList, epsilon, minpts, calDistance=calLatLngDistance):
    ind2directDenseReachInde = {}
    for indF in xrange(len(latlngList)):
        lngF , latF = latlngList[indF][0:2]
        for indT in xrange(len(latlngList)):
            lngT, latT = latlngList[indT][0:2]
            #优化 ，防止重复计算，
            if latT in ind2directDenseReachInde and latF in ind2directDenseReachInde[latT]:
                ind2directDenseReachInde[indF].add(indT)
                continue
            if latT in ind2directDenseReachInde and latF not in ind2directDenseReachInde[latT]:
                continue
            dis = calDistance(lngF, latF, lngT, latT)
            if dis <= epsilon:
                if indF not in ind2directDenseReachInde:
                    ind2directDenseReachInde[indF] = set()
                ind2directDenseReachInde[indF].add(indT)
    # 已添加进簇的 ind
    clustered = set()

    # 返回所有簇
    clusters = []
    for ind in ind2directDenseReachInde:
        if len(ind2directDenseReachInde[ind]) >= minpts and ind not in clustered:
            # ind 为核心点 。遍历该核心点邻域内的所有核心点,开始expand形成新簇
            clustered.add(ind)

            # 一个新簇
            newCluster = []
            newCluster.append(ind)
            # 扩充队列
            tmpClusterQueue = []
            tmpClusterQueue.append(ind)
            while len(tmpClusterQueue) != 0:
                curInd = tmpClusterQueue.pop()
                for inInd in ind2directDenseReachInde[curInd]:
                    if len(ind2directDenseReachInde[inInd]) < minpts and inInd not in clustered:
                        clustered.add(inInd)
                        newCluster.append(inInd)
                        # inInd 为边界点，无法扩充
                    elif len(ind2directDenseReachInde[inInd]) >= minpts and inInd not in clustered:
                        clustered.add(inInd)
                        newCluster.append(inInd)
                        # inInd 为核心点，开始扩充 ,加入tmpClusterQueue队列
                        tmpClusterQueue.append(inInd)
            if len(newCluster) != 0:
                clusters.append(newCluster)

    noises = [key for key in ind2directDenseReachInde.keys() if key not in clustered]
    clustersP = [[latlngList[ind] for ind in oneCind] for oneCind in clusters]
    noisesP = [latlngList[ind] for ind in noises]
    return {"clusters": clustersP, "noises": noisesP}


# 返回簇中心经纬度,簇内点个数,半径，ssid分布
# ret: {"pointNums": 18, "ssid": "Lucifer", "ssidSta": {"Lucifer": 18}, "centerLng": 22.75629348625081, "radius": 52, "centerLat": 113.87071855056054}
def getWifiCenterLatLng(clusteredlatlngList, useRssi=True, calDistance=calLatLngDistance):
    weightFunc = lambda rssi: 1
    if useRssi:
        weightFunc = lambda rssi: 1 / (1 + math.exp(-1 * rssi))
    weigthLats = [one[1] * weightFunc(rssi=one[2]) for one in clusteredlatlngList  ]
    weigthLngs = [one[0] * weightFunc(rssi=one[2]) for one in clusteredlatlngList  ]
    weigths = [weightFunc(rssi=one[2]) for one in clusteredlatlngList ]
    ssidMap = {}
    for ssid in [one[3] for one in clusteredlatlngList if one[2] <= 0.0]:
        if ssid in ssidMap:
            ssidMap[ssid] = ssidMap[ssid] + 1
        else:
            ssidMap[ssid] = 1
    ssidSort = sorted(ssidMap.items(), cmp=lambda x, y: cmp(x[1], y[1]), reverse=True)
    # 选择簇内离中心店最远点的距离作为radius
    centerLat = sum(weigthLats) / sum(weigths)
    centerLng = sum(weigthLngs) / sum(weigths)
    radius = max(
        [int(calDistance(centerLng , centerLat, one[0], one[1])) for one in clusteredlatlngList  ])
    return {"centerLat": centerLat, "centerLng": centerLng, "pointNums": len(weigths), "radius": radius, "ssidSta": ssidMap ,"ssid": ssidSort[0][0]}

import time
if __name__=="__main__":
    # data = '["26.887899,100.22496,-86,@往返免费WIFI","24.963003,102.795415,-92,wjscy","24.963316,102.795008,-80,wjscy","24.962798,102.795044,-89,wjscy","24.962619,102.795353,-80,wjscy","24.962818,102.795803,-79,wjscy","24.962562,102.79547,-80,wjscy","24.962182,102.796018,-87,wjscy","24.962317,102.793327,-88,wjscy","24.962558,102.795479,-79,wjscy","24.962459,102.79365,-81,wjscy","24.96322,102.79493,-71,wjscy","24.963074,102.795211,-81,wjscy","24.962837,102.795955,-90,wjscy","24.962665,102.795227,-81,wjscy","24.962592,102.795009,-77,wjscy","24.962726,102.794972,-94,wjscy","24.962936,102.795237,-81,wjscy","24.962153,102.796074,-88,wjscy","24.962583,102.796847,-82,wjscy","24.962651,102.795084,-89,wjscy","22.919362,114.079854,-76,wjscy","22.922532,114.081617,0,wjscy","22.919008,114.079683,-79,wjscy","22.919048,114.079651,-70,wjscy","22.919363,114.079819,-77,wjscy","22.919036,114.079666,-72,wjscy","22.919082,114.079675,-52,wgmzl"]'
    # data = '["23.09985,113.27106,-91,YINGJIA-MAN","23.11680985854638,113.27764549750246,-89,dd-wrt","23.116883,113.275798,-81,dd-wrt","23.09985,113.27111,-82,YINGJIA-MAN","23.116138,113.276444,-87,dd-wrt","23.099870210208056,113.2710438138049,-83,YINGJIA-MAN","23.099748,113.270952,-92,YINGJIA-MAN","23.115945172608818,113.27723805092188,-88,dd-wrt","23.115995250060173,113.27730812068559,-75,dd-wrt","23.09985,113.27106,-80,YINGJIA-MAN","23.09985,113.27106,-88,YINGJIA-MAN","23.09983,113.27111,-91,YINGJIA-MAN","0.0,0.0,-87,dd-wrt","23.09988,113.27107,-86,YINGJIA-MAN","23.0999,113.27107,-86,YINGJIA-MAN","23.115842262239404,113.27739919700821,-84,dd-wrt","23.09985,113.27111,-84,YINGJIA-MAN","23.09985,113.27111,-86,YINGJIA-MAN","23.11607,113.27857,-90,dd-wrt","23.09985,113.27106,-82,YINGJIA-MAN","23.116532,113.276853,-95,dd-wrt","23.099755,113.270947,-90,YINGJIA-MAN","23.116330347254888,113.27726009761808,-72,dd-wrt","23.09983,113.27111,-85,YINGJIA-MAN","23.116475,113.276866,-91,dd-wrt","23.11599,113.27731,-75,dd-wrt","23.09979318132812,113.27104881223444,-82,YINGJIA-MAN","23.115995250060173,113.27730812068559,-79,dd-wrt","23.117705,113.27821,-91,dd-wrt","23.099846,113.271094,-85,YINGJIA-MAN","23.09990524018388,113.27106983483867,-75,YINGJIA-MAN","23.116065991160262,113.2769527904958,-87,dd-wrt","23.09988723094779,113.27106683147439,-77,YINGJIA-MAN","23.116064313306968,113.27735116619803,-83,dd-wrt","23.11596020234892,113.27726707941899,-80,dd-wrt","23.09979,113.27095,-67,YINGJIA-MAN","23.145552,113.252116,-68,POGXEV","23.145564,113.252122,-79,POGXEV","23.099757,113.270946,-86,YINGJIA-MAN","23.09987,113.27105,-81,YINGJIA-MAN","23.099857246006383,113.27111286258284,-79,YINGJIA-MAN","23.09985,113.27111,-87,YINGJIA-MAN","23.09985,113.27111,-79,YINGJIA-MAN"]'
    # data = '["34.093195,119.271926,-59,@往返免费WIFI","42.257836,118.959471,-87,@往返免费WIFI","29.644792,119.536459,-90,@往返免费WIFI","34.093195,119.271926,-63,@往返免费WIFI","22.632214,114.796985,-87,@往返免费WIFI","42.656129,118.928813,-86,@往返免费WIFI","35.424636,119.570379,-89,@往返免费WIFI","42.258581,118.959123,-89,@往返免费WIFI","22.677175,114.756225,-82,@往返免费WIFI","42.2581365097986,118.95937380697049,-77,@往返免费WIFI","42.25807536784416,118.95930659371561,-76,@往返免费WIFI","42.25809456284507,118.95939687141887,-77,@往返免费WIFI","42.934974,118.998112,-73,@往返免费WIFI","42.258468,118.959108,-89,@往返免费WIFI","37.515736,105.218272,-85,@往返免费WIFI","42.257807010900116,118.95959343139984,-77,@往返免费WIFI","35.45435,119.581728,-73,@往返免费WIFI","43.127753,119.058355,-84,@往返免费WIFI","22.67717,114.756226,-92,@往返免费WIFI","29.566112,103.448917,-91,@往返免费WIFI","42.260613,118.957556,-76,@往返免费WIFI","29.566071,103.446624,-80,@往返免费WIFI","43.126048,119.061094,-92,@往返免费WIFI","35.466103,119.537081,-70,@往返免费WIFI","22.635272,114.789278,-83,@往返免费WIFI","40.001749,119.75793,-90,@往返免费WIFI","29.580378,103.425457,-88,@往返免费WIFI","29.56638,103.450166,-77,@往返免费WIFI","22.662188,114.744941,-91,@往返免费WIFI","35.455605,119.583898,-54,@往返免费WIFI","30.854259,104.404666,-82,@往返免费WIFI","35.413564,119.565336,-75,@往返免费WIFI","37.839211,112.58376,-75,@往返免费WIFI","42.273487,118.91238,-85,@往返免费WIFI","35.412441,119.563376,-84,@往返免费WIFI","30.31094983698648,120.15446106282205,-73,@往返免费WIFI","22.674586,114.754897,-65,@往返免费WIFI","29.567807,103.454625,-88,@往返免费WIFI","35.408863,119.549128,-71,@往返免费WIFI","35.388568,119.548782,-78,@往返免费WIFI","35.45246,119.579977,-81,@往返免费WIFI","29.575752,103.407748,-54,@往返免费WIFI","29.566544,103.450106,-72,@往返免费WIFI","42.258471,118.959222,-77,@往返免费WIFI","35.454436,119.581635,-69,@往返免费WIFI","37.839557,112.586168,-85,@往返免费WIFI","37.51512928436265,105.21861911357462,-83,@往返免费WIFI","42.244787,118.927093,-90,@往返免费WIFI","37.86016,112.576778,-71,@往返免费WIFI","22.675142800745473,114.75416180389378,-70,@往返免费WIFI","35.464094,119.55308,-47,@往返免费WIFI","43.125894,119.061218,-90,@往返免费WIFI","42.25781,118.95959,-76,@往返免费WIFI","42.2583,118.95873,-77,@往返免费WIFI","35.468509,119.592429,-85,@往返免费WIFI","33.519617,119.156243,-64,@往返免费WIFI","29.583027,103.408764,-92,@往返免费WIFI","35.468542,119.592449,-88,@往返免费WIFI","30.31094,120.15446,-73,@往返免费WIFI","42.25810258631146,118.95940790618047,-77,@往返免费WIFI","35.413589,119.565343,-77,@往返免费WIFI","33.876959,119.230656,-55,@往返免费WIFI","35.413663,119.565137,-78,@往返免费WIFI","37.857008,112.58238,-86,@往返免费WIFI","40.004382,119.765057,-76,@往返免费WIFI","42.264371,118.952112,-81,@往返免费WIFI","30.196778,120.134112,-76,@往返免费WIFI","35.413563,119.56531,-64,@往返免费WIFI","35.413576,119.565298,-67,@往返免费WIFI","29.575527,103.40779,-80,@往返免费WIFI","42.263253,118.96725,-90,@往返免费WIFI","42.935025,118.998085,-68,@往返免费WIFI","29.581915,103.4307,-80,@往返免费WIFI","35.412233,119.563275,-82,@往返免费WIFI","40.00579,119.766852,-95,@往返免费WIFI","40.05749,119.86982,-87,@往返免费WIFI","42.257899,118.959354,-78,@往返免费WIFI"]'
    # data = '["21.880647,111.932895,-88,TP-LINK_BBF9","21.880524,111.932948,-93,TP-LINK_BBF9","21.880659,111.932905,-80,TP-LINK_BBF9","21.880659,111.932905,-88,TP-LINK_BBF9","21.880705,111.932928,-74,TP-LINK_BBF9","21.880705,111.932928,-73,TP-LINK_BBF9","21.880707,111.93293,-74,TP-LINK_BBF9","21.880705,111.932928,-70,TP-LINK_BBF9","21.880707,111.93293,-73,TP-LINK_BBF9","21.880705,111.932928,-72,TP-LINK_BBF9","21.880702,111.932933,-73,TP-LINK_BBF9","21.880707,111.93293,-72,TP-LINK_BBF9"]'
    # data = '["39.90817,116.42764,-76,CU_7K0B","39.907951,116.427759,-86,CU_7K0B","39.908057,116.427698,-86,CU_7K0B","39.907602098215236,116.42757921120442,-61,CU_7K0B","39.907586,116.427772,-70,CU_7K0B","39.907549,116.427854,-67,CU_7K0B","39.90767,116.427709,-66,CU_7K0B","39.907538,116.427719,-70,CU_7K0B","39.908059,116.427616,-68,CU_7K0B","39.908057,116.42758,-60,CU_7K0B","39.908071,116.427648,-68,CU_7K0B"]'
    # dataTransform = [ [ item if i==3  else float(item) for i,item in enumerate(one.split(",")) ] for one in json.loads(data) ]
    # clusters = dbscan(dataTransform,200,5,calDistance=calLatLngDistance)
    # print clusters
    # for c in clusters:
    #     if c == "clusters":
    #         for oneC in clusters[c]:
    #             print oneC
    #             print "\t",getWifiCenterLatLng(oneC, True)
    # exit(0)
    # 模拟输入
    tstart= time.time()
    for line in open("/Users/leichen/Downloads/wifi_sample_more.txt"):
        segs = line.strip('\n').split('\t')
        collect_latlngrssi = segs[0]
        epsilon = float(segs[1])
        minpts = int(segs[2])
        otherWithData = "\t".join(segs[3:])
        # try:
        wifiLatLngs = [[item if i >= 3 else float(item) for i, item in enumerate(one.split(","))] for one in json.loads(collect_latlngrssi)][0:500]
        clusterRet = dbscan(wifiLatLngs, epsilon, minpts, calDistance=calLatLngDistance)
        if len(clusterRet['clusters']) == 0:
            # 防止点数过少、等原因导致无聚类结果 。 epsilon=200 minpts=1 保底
            clusterRet = dbscan(wifiLatLngs, 200, 1, calDistance=calLatLngDistance)
        sortedClusterCenter = sorted([getWifiCenterLatLng(oneC) for oneC in clusterRet['clusters']],
                                     cmp=lambda x, y: cmp(x["pointNums"], y["pointNums"]), reverse=True)
        print json.dumps(sortedClusterCenter[0]) \
                  + "\t" + json.dumps(sortedClusterCenter) \
                  + "\t" + otherWithData
        # except Exception,e:
        #     print collect_latlngrssi
        #     print e
        # exmaple : 24.9632179058,102.794930838,15,wjscy:15	[[24.963217905806083, 102.79493083844251, 15, "wjscy:15"], [22.9190819999994, 114.07967499999961, 6, "wjscy:5|wgmzl:1"]]	{"clusters": [[[24.963003, 102.795415, -92.0, "wjscy"], [24.962798, 102.795044, -89.0, "wjscy"], [24.962619, 102.795353, -80.0, "wjscy"], [24.962818, 102.795803, -79.0, "wjscy"], [24.962562, 102.79547, -80.0, "wjscy"], [24.962558, 102.795479, -79.0, "wjscy"], [24.963074, 102.795211, -81.0, "wjscy"], [24.962665, 102.795227, -81.0, "wjscy"], [24.962936, 102.795237, -81.0, "wjscy"], [24.963316, 102.795008, -80.0, "wjscy"], [24.96322, 102.79493, -71.0, "wjscy"], [24.962592, 102.795009, -77.0, "wjscy"], [24.962726, 102.794972, -94.0, "wjscy"], [24.962651, 102.795084, -89.0, "wjscy"], [24.962837, 102.795955, -90.0, "wjscy"]], [[22.919362, 114.079854, -76.0, "wjscy"], [22.919008, 114.079683, -79.0, "wjscy"], [22.919048, 114.079651, -70.0, "wjscy"], [22.919363, 114.079819, -77.0, "wjscy"], [22.919036, 114.079666, -72.0, "wjscy"], [22.919082, 114.079675, -52.0, "wgmzl"]]], "noises": [[24.962182, 102.796018, -87.0, "wjscy"], [24.962317, 102.793327, -88.0, "wjscy"], [24.962459, 102.79365, -81.0, "wjscy"], [24.962153, 102.796074, -88.0, "wjscy"], [24.962583, 102.796847, -82.0, "wjscy"], [22.922532, 114.081617, 0.0, "wjscy"]]}
    tend = time.time()
    print tend-tstart


